Geography of Covid Vaccine Hesitancy in the U.S.

Vaccine Line in Boston Massachusets, March 2021 # Introduction Back in March, I created my PUG Shiny project with a focus on COVID vaccine distribution. At the time, the issue facing the united states was primarily one of distribution: people who wanted the vaccine couldn’t get it. Eligibility was restricted state by state to allocate a precious resource. Even people who were eligible needed to get lucky and navigate online sign up forms, vouching for fewer spots than people. Here’s my Shiny App, improved and updated with the latest data:

`

###############
# import data #
###############
both <- read_csv("both.csv")
## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   jurisdiction = col_character(),
##   week = col_character(),
##   first_dose_m = col_double(),
##   second_dose_m = col_double(),
##   all_doses_m = col_double(),
##   cumulative_doses_m = col_double(),
##   first_dose_p = col_double(),
##   second_dose_p = col_double(),
##   all_doses_p = col_double(),
##   cumulative_doses_p = col_double(),
##   first_dose_both = col_double(),
##   second_dose_both = col_double(),
##   all_doses_both = col_double(),
##   cumulative_doses_both = col_double()
## )
both[, "week"]<-as.Date(both$week, "%m/%d/%Y")
both <- both[order(both$"week"),] 


###################################################
# CHOICE VALUES AND LABELS#
###################################################
vacc_choice_values <- c("_m","_p","_both")
vacc_choice_names <- c("Moderna","Pfizer","Both Moderna and Pfizer")
names(vacc_choice_values) <- vacc_choice_names

dose_choice_values <- c("all_doses", "cumulative_doses")
dose_choice_names <- c("Total New Doses this Week", "Cumulative Doses to Date")
names(dose_choice_values) <- dose_choice_names

state_choice_values <- unique(both$jurisdiction)
state_choice_names <- sapply(state_choice_values, FUN = str_to_title)
names(state_choice_values) <- state_choice_names

week_choices <- unique(both$week)






############
#    ui    #
############
ui <- fluidPage(
  
  titlePanel(strong("Weekly Vaccine Allocations by State")),
  
  tabPanel(
    h1(style = "font-family:Impact",
       "Map"),
    sidebarLayout(
      sidebarPanel(
        #SELECT VACCINE BRAND
        selectInput(inputId = "vaccvar", 
                    label = "Choose a vaccine type to plot:", 
                    choices = vacc_choice_values, 
                    selected = "_p"),
        #SELECT DOSE NUMBER
        selectInput(inputId = "dosevar", 
                    label = "Choose a dose number to plot:", 
                    choices = dose_choice_values, 
                    selected = "first_"),
        #SELECT MULTIPLE STATES
        selectizeInput(inputId = "state_name", 
                       label = "Pick states to visualize:", 
                       choices = state_choice_values, 
                       selected = "california", 
                       multiple = TRUE),
        #SLIDER DATE
        sliderTextInput(inputId = "slider", 
                        label = h4(tags$b("Choose a week to plot:")), 
                        choices = week_choices, 
                        selected = "2020-12-21", 
                        grid = TRUE)),
      mainPanel(plotOutput(outputId = "distPlot"))
    )
  )
)
#~~~~~~~~~~~#
# END OF UI #
#~~~~~~~~~~~#







############
# server   #
############

server <- function(input, output) {
  
  output$distPlot <- renderPlot({
    
    dat_vacc <- both%>%
      select(c("jurisdiction", "week", contains(input$vaccvar)))%>%
      select(c("jurisdiction", "week", contains(input$dosevar)))%>%
      filter(week == input$slider)
    
    p <-ggplot(data = dat_vacc %>%filter(jurisdiction %in% input$state_name), aes(map_id = jurisdiction)) +    
      geom_map(aes(fill = get(paste(input$dosevar, input$vaccvar, sep = ""))), map = fifty_states) + 
      
      ###WANT A MAX VALUE FOR FIRST, SECOND, ALL, AND CUMULATIVE
      
      #scale_fill_gradient(low="white", high="dark green") +
      #scale_fill_gradient(low="blue", mid = "white", high="red") +
      scale_fill_gradientn(colors = rainbow(3)) +
      scale_fill_gradientn(colors = c("white", "blue", "red")) +
      #scale_fill_continuous(colors = c(rainbow(5))) +
      expand_limits(x = fifty_states$long, y = fifty_states$lat) +
      
      coord_map() +
      scale_x_continuous(breaks = NULL) + 
      scale_y_continuous(breaks = NULL) +
      
      labs(x = "", y = "", 
           title = str_wrap(
             paste(
               names(which(vacc_choice_values == input$vaccvar)), 
               names(which(dose_choice_values == input$dosevar)),
               " to", 
               paste(str_to_title(unlist(input$state_name)), collapse=', ')
             ), 
             80),
           subtitle = paste("Week of ", input$slider))+
      
      theme(legend.position = "left", 
            legend.direction = "vertical",
            panel.background = element_rect(fill = alpha("light blue", 0.3)), 
            legend.key.size = unit(1.5, "cm"),
            legend.key.width = unit(0.5,"cm"),
            text=element_text(size=16,  family="serrif"), 
            plot.title = element_text(hjust = 0.5),
            plot.subtitle = element_text(hjust = 0.5))+
      
      labs(fill = "Doses")
    
    
    plot_scale <- function(dosevar) {
      
      if(input$dosevar== "combined_doses" & input$vaccvar == "_both")
        scale_fill_gradient(low="blue", high="red", limits = c(0, 40000000))
      
      else if(input$vaccvar == "_both" & input$dosevar != "cumulative_doses")
        scale_fill_gradient(low="blue", high="red", limits = c(0, 2000000))
      
      else if(input$dosevar%in% c("first_dose", "second_dose"))
        scale_fill_gradient(low="blue", high="red", limits = c(0, 650000))
      
      else if(input$dosevar == "all_doses")
        scale_fill_gradient(low="blue", high="red", limits = c(0, 750000))
      
      else if(input$dosevar== "cumulative_doses" & input$vaccvar == "_both")
        scale_fill_gradient(low="blue", high="red", limits = c(0, 22000000))
      else if(input$dosevar== "cumulative_doses")
        scale_fill_gradient(low="blue", high="red", limits = c(0, 22000000))
    }
    
    p + plot_scale()
  })
  
}

shinyApp(ui = ui, server = server)
## 
## Listening on http://127.0.0.1:3973
## Scale for 'fill' is already present. Adding another scale for 'fill', which
## will replace the existing scale.
## Scale for 'fill' is already present. Adding another scale for 'fill', which
## will replace the existing scale.

Motivation

Over the past six weeks, the supply chain has greatly improved, and access is no longer such a pressing issue. The issues facing the country have evolved: the issue now is acheiving heard immunity, and getting everyone vaccinated. In March, I didn’t know anyone who had gotten a vaccine. But now in April, I don’t know anyone who hasn’t gotten one. I wanted to map out where in the country people are the most hesitant to get vaccinated. Knowing the most resistant geographies could be useful information for government organizations to focus their vaccine advocacy efforts.

Anti-Vax Protesters in Texas

Mapping County-Level Vaccine Hesitancy

I investigated a CDC published dataset titled “Vaccine Hesitancy for COVID-19: County and local estimates.” The dataset was produced by the ASPE, the Office of the Assistant Secretary for Planning and Evaluation, of the U.S. Department of Health and Human Services. Specifically, the data was collected via the U.S. Census Bureau’s Household Pulse Survey.

H2: Variables of Interest

The dataset contained several demographic variables measured for each U.S. county. The variables of interest were “Estimated Hesitant” and “Estimated Strongly Hesitant”, which measured the percent of adults in each county who indicated either hesitant, or strongly hesitant, respectively, to get vaccinated. I assume that survey respondants are all adults who have not yet received any covid vaccine.

H2: Data Wrangling

H3: Fauci with a headache, or …

fauchi

H3: Data Wrangling

{r} # my_data0 <- my_data %>% # select(!c("Percent non-Hispanic Native Hawaiian/Pacific Islander", "FIPS Code"))%>% # rename(region = State, subregion = `County Name`) %>% # mutate(region = tolower(region), subregion = tolower(subregion))%>% # rename(est_hesitant = "Estimated hesitant")%>% # rename(est_strong_hesitant = "Estimated strongly hesitant")%>% # rename(svi_cat = "SVI Category" )%>% # rename(svi = "Social Vulnerability Index (SVI)" )%>% # rename(pct_full_vaxed = "Percent adults fully vaccinated against COVID-19 as of 3/30/2021") # # my_data0$subregion <- str_remove(my_data0$subregion, " county.*") # my_data0$subregion <- str_remove(my_data0$subregion, ",.*") # # # census0 <- census%>% # mutate(subregion = gsub("\\.", "", subregion))%>% # mutate(subregion = tolower(subregion))%>% # mutate(region = gsub(".*, ", "",subregion)) %>% # mutate(subregion = gsub(", .*", "", subregion))%>% # rename(pop_est = "2019_pop_est") # census0$subregion <- str_remove(census0$subregion, " county.*") # # my_data0 <- my_data0 %>% # inner_join(census0, by =c("subregion", "region")) # head(my_data0, 100) #

H2: Making the maps

H3: Code

Creating mapable data set

{r} # usa_counties <- map_data(map = "county", region = ".") # # my_data0_map <- my_data0 %>% # inner_join(usa_counties, by =c("subregion", "region")) #

Mapping of hesitancy

{r,eval= FALSE} # ggplot(my_data0_map, aes(x = long, y = lat, group = group, fill = est_hesitant)) + # geom_polygon(color = "white", size = 0.05) + # theme_void() + # coord_fixed(ratio = 1.3) + # labs(fill = "Proportion of residents hesitant to be vaccinated") + # theme(legend.position="bottom")+ # scale_fill_distiller(palette = "Spectral") #

Mapping of strong hesitancy

{r, eval=FALSE} # ggplot(my_data0_map, aes(x = long, y = lat, group = group, fill = est_strong_hesitant)) + # geom_polygon(color = "white", size = 0.05) + # theme_void() + # coord_fixed(ratio = 1.3) + # labs(fill = "Proportion of residents STRONGLY hesitant to be vaccinated") + # theme(legend.position="bottom")+ # scale_fill_distiller(palette = "Spectral") #

H3: Hesitancy Map

{r,echo= FALSE} # ggplot(my_data0_map, aes(x = long, y = lat, group = group, fill = est_hesitant)) + # geom_polygon(color = "white", size = 0.05) + # theme_void() + # coord_fixed(ratio = 1.3) + # labs(fill = "Proportion of residents hesitant to be vaccinated") + # theme(legend.position="bottom")+ # scale_fill_distiller(palette = "Spectral") #

H3: Strong Hesitancy Map

{r, echo= FALSE} # ggplot(my_data0_map, aes(x = long, y = lat, group = group, fill = est_strong_hesitant)) + # geom_polygon(color = "white", size = 0.05) + # theme_void() + # coord_fixed(ratio = 1.3) + # labs(fill = "Proportion of residents STRONGLY hesitant to be vaccinated") + # theme(legend.position="bottom")+ # scale_fill_distiller(palette = "Spectral") #

H3: zooming in on distinct areas

{r midwest} # ggplot() + # geom_polygon(data = my_data0_map, aes(x = long, y = lat, group = group, fill = est_strong_hesitant), color = "white", size = 0.05) + # geom_polygon(data = usa_states, aes(x = long, y = lat, group = group),fill="NA", colour = "black") + # theme_void() + # coord_fixed(xlim = c(-110, -90), ylim = c(42, 50), ratio = 1.3) + # labs(fill = "Proportion of residents STRONGLY hesitant to be vaccinated") + # theme(legend.position="bottom")+ # scale_fill_distiller(palette = "Spectral") # # remove legend! # guides(fill = FALSE) #

{r tenn} # ggplot() + # geom_polygon(data = my_data0_map, aes(x = long, y = lat, group = group, fill = est_hesitant), color = "white", size = 0.05) + # geom_polygon(data = usa_states, aes(x = long, y = lat, group = group),fill="NA", colour = "black") + # theme_void() + # coord_fixed(xlim = c(-95, -75), ylim = c(30, 40), ratio = 1.3) + # labs(fill = "Proportion of residents STRONGLY hesitant to be vaccinated") + # theme(legend.position="bottom")+ # scale_fill_distiller(palette = "Spectral") # # remove legend! # guides(fill = FALSE) #

H1: Analysis of Maps

I observe a striking visual pattern in the mapping of county hesitancy levels. In many places on the map, we see areas that dont

IDK yet but theres something about how state lines really impact fill. I ran a cluser analysis to see if I was making this trend up or not

Try out clustering # {r, warning=FALSE, message=FALSE} # my_data0_map_cluster <- my_data0_map %>% # select(est_hesitant, long, lat) # # # set.seed(15) # library(mclust) # county_clusts <- my_data0_map_cluster %>% # kmeans(centers = 48)%>% # fitted("classes")%>% # as.character() # # my_data0_map_cluster <- my_data0_map_cluster %>% mutate(cluster = county_clusts) # # my_data0_map_cluster %>% ggplot(aes(x = long, y = lat)) + # geom_point(aes(color = cluster), alpha = 0.5)+ # coord_fixed(ratio = 1.3) #

{r, echo=FALSE} # #usa_states <- map_data(map = "state", region = ".") # # ggplot()+ # geom_point(data = my_data0_map_cluster, aes(x = long, y = lat, color = cluster), alpha = 0.5)+ # geom_polygon(data = usa_states, aes(x = long, y = lat, group = group),fill="NA", colour = "black") + # theme_void()+ # theme(legend.position="none") #

How to make tabs

Bulleted list

You can make a bulleted list like this:

  • item 1
  • item 2
  • item 3

Numbered list

You can make a numbered list like this

  1. First thing I want to say
  2. Second thing I want to say
  3. Third thing I want to say